{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Content under Creative Commons Attribution license CC-BY 4.0, code under BSD 3-Clause License © 2018 D. Koehn, notebook style sheet by L.A. Barba, N.C. Clementi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<link href=\"https://fonts.googleapis.com/css?family=Merriweather:300,300i,400,400i,700,700i,900,900i\" rel='stylesheet' >\n",
       "<link href=\"https://fonts.googleapis.com/css?family=Source+Sans+Pro:300,300i,400,400i,700,700i\" rel='stylesheet' >\n",
       "<link href='http://fonts.googleapis.com/css?family=Source+Code+Pro:300,400' rel='stylesheet' >\n",
       "<style>\n",
       "\n",
       "@font-face {\n",
       "    font-family: \"Computer Modern\";\n",
       "    src: url('http://mirrors.ctan.org/fonts/cm-unicode/fonts/otf/cmunss.otf');\n",
       "}\n",
       "\n",
       "\n",
       "#notebook_panel { /* main background */\n",
       "    background: rgb(245,245,245);\n",
       "}\n",
       "\n",
       "div.cell { /* set cell width */\n",
       "    width: 800px;\n",
       "}\n",
       "\n",
       "div #notebook { /* centre the content */\n",
       "    background: #fff; /* white background for content */\n",
       "    width: 1000px;\n",
       "    margin: auto;\n",
       "    padding-left: 0em;\n",
       "}\n",
       "\n",
       "#notebook li { /* More space between bullet points */\n",
       "margin-top:0.5em;\n",
       "}\n",
       "\n",
       "/* draw border around running cells */\n",
       "div.cell.border-box-sizing.code_cell.running { \n",
       "    border: 1px solid #111;\n",
       "}\n",
       "\n",
       "/* Put a solid color box around each cell and its output, visually linking them*/\n",
       "div.cell.code_cell {\n",
       "    background-color: rgb(256,256,256); \n",
       "    border-radius: 0px; \n",
       "    padding: 0.5em;\n",
       "    margin-left:1em;\n",
       "    margin-top: 1em;\n",
       "}\n",
       "\n",
       "\n",
       "div.text_cell_render{\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    line-height: 140%;\n",
       "    font-size: 110%;\n",
       "    width:680px;\n",
       "    margin-left:auto;\n",
       "    margin-right:auto;\n",
       "}\n",
       "\n",
       "/* Formatting for header cells */\n",
       ".text_cell_render h1 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-style:regular;\n",
       "    font-weight: bold;    \n",
       "    font-size: 250%;\n",
       "    line-height: 100%;\n",
       "    color: #004065;\n",
       "    margin-bottom: 1em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       ".text_cell_render h2 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: bold; \n",
       "    font-size: 180%;\n",
       "    line-height: 100%;\n",
       "    color: #0096d6;\n",
       "    margin-bottom: 0.5em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       "\n",
       ".text_cell_render h3 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "\tfont-size: 150%;\n",
       "    margin-top:12px;\n",
       "    margin-bottom: 3px;\n",
       "    font-style: regular;\n",
       "    color: #008367;\n",
       "}\n",
       "\n",
       ".text_cell_render h4 {    /*Use this for captions*/\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: 300; \n",
       "    font-size: 100%;\n",
       "    line-height: 120%;\n",
       "    text-align: left;\n",
       "    width:500px;\n",
       "    margin-top: 1em;\n",
       "    margin-bottom: 2em;\n",
       "    margin-left: 80pt;\n",
       "    font-style: regular;\n",
       "}\n",
       "\n",
       ".text_cell_render h5 {  /*Use this for small titles*/\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    font-weight: regular;\n",
       "    font-size: 130%;\n",
       "    color: #e31937;\n",
       "    font-style: italic;\n",
       "    margin-bottom: .5em;\n",
       "    margin-top: 1em;\n",
       "    display: block;\n",
       "}\n",
       "\n",
       ".text_cell_render h6 { /*use this for copyright note*/\n",
       "    font-family: 'Source Code Pro', sans-serif;\n",
       "    font-weight: 300;\n",
       "    font-size: 9pt;\n",
       "    line-height: 100%;\n",
       "    color: grey;\n",
       "    margin-bottom: 1px;\n",
       "    margin-top: 1px;\n",
       "}\n",
       "\n",
       "    .CodeMirror{\n",
       "            font-family: \"Source Code Pro\";\n",
       "\t\t\tfont-size: 90%;\n",
       "    }\n",
       "/*    .prompt{\n",
       "        display: None;\n",
       "    }*/\n",
       "\t\n",
       "    \n",
       "    .warning{\n",
       "        color: rgb( 240, 20, 20 )\n",
       "        }  \n",
       "</style>\n",
       "<script>\n",
       "    MathJax.Hub.Config({\n",
       "                        TeX: {\n",
       "                           extensions: [\"AMSmath.js\"], \n",
       "                           equationNumbers: { autoNumber: \"AMS\", useLabelIds: true}\n",
       "                           },\n",
       "                tex2jax: {\n",
       "                    inlineMath: [ ['$','$'], [\"\\\\(\",\"\\\\)\"] ],\n",
       "                    displayMath: [ ['$$','$$'], [\"\\\\[\",\"\\\\]\"] ]\n",
       "                },\n",
       "                displayAlign: 'center', // Change this to 'center' to center equations.\n",
       "                \"HTML-CSS\": {\n",
       "                    styles: {'.MathJax_Display': {\"margin\": 4}}\n",
       "                }\n",
       "        });\n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Execute this cell to load the notebook's style sheet, then ignore it\n",
    "from IPython.core.display import HTML\n",
    "css_file = '../../style/custom.css'\n",
    "HTML(open(css_file, \"r\").read())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Equidistant Cartesian grid\n",
    "\n",
    "In [lecture 1](https://github.com/daniel-koehn/Theory-of-seismic-waves-II/tree/master/01_Analytical_solutions) we reviewed the forward problems to describe wave propagation in isotropic elastic and acoustic media. Furthermore, analytical solutions for homogeneous 3D, 2D and 1D acoustic media are derived to verify the accuarcy of the finite difference codes we develop later in the course.\n",
    "\n",
    "Before introducing the basics of finite-difference modelling, I want to give an overview of some concepts to generate discrete subsurface models."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## From the continous to discrete world  \n",
    "\n",
    "The partial differential equations, which we developed before describe wave propagation in a continous medium. However, if we want to solve these equations in a computer, we have to discretize the continuum, due to the limitation that a computer can only compute discrete problems. The question is how fine do we have to discretize the continuum in order to get accurately seismic modelling results?\n",
    "\n",
    "Do we have to discretize down to the molecular scale? For example 12 grams of the carbon isotop C-12 contains $6.022140857x10^{23}$ atoms. Let's assume that we can represent the coordinates of the atoms in each direction by a double precision floating-point number which requires 8 Byte, so in totally we need ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$1.4x10^{13}$ TB would lead to some serious memory problems when using such a fine mesh, because the [NEC-cluster](https://www.rz.uni-kiel.de/de/angebote/hiperf/nec-linux-cluster) at the computing centre at Christian-Albrechts-University Kiel only has a total of 4 TB of RAM. We have not even considered the required runtime for such a seismic modelling run.\n",
    "\n",
    "Fortunately, we do not have to discretize the universe down to the atomic or moleculare scale in order to get reasonable modelling results. In a later lecture, we will derive criteria to find the best compromise between accuracy of the modelling results and level of model discretization."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Discretization of a 2D model\n",
    "\n",
    "Let's start with the 2D acoustic wave equation:\n",
    "\n",
    "\\begin{equation}\n",
    "\\frac{\\partial^2 P(x,z,t)}{\\partial t^2} - Vp^2(x,z) \\biggl(\\frac{\\partial^2 P(x,z,t)}{\\partial x^2} + \\frac{\\partial^2 P(x,z,t)}{\\partial z^2} \\biggl)= f \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "In a first step to solve this problem using finite-differences we have to discretize the pressure wavefield $P(x,z,t)$ and $V_{p0}(x,z)$ at a specific time t on a **equidistant 2D Cartesian spatial grid**:\n",
    "\n",
    "<img src=\"../images/2D-grid_cart_ac.png\" width=\"75%\">\n",
    "\n",
    "where we define the positions (x,z) of the discrete points in the model by\n",
    "\n",
    "\\begin{align}\n",
    "x &= i*dh,\\nonumber \\\\\n",
    "z &= j*dh,\\nonumber \\\\\n",
    "\\end{align}\n",
    "\n",
    "where $dh$ is the spatial grid point distance. We denote the value of the pressure wavefield and P-wave velocity at a discrete point by $P_{j,i}$ and $Vp_{j,i}$, respectively. These are 2D matrices, so we can simply store them in 2D arrays. Let's try to do this in Python."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "code_folding": [
     0
    ]
   },
   "outputs": [],
   "source": [
    "# Import Libraries\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Here, I introduce a new library, which is useful \n",
    "# to define the fonts and size of a figure in a notebook\n",
    "from pylab import rcParams\n",
    "\n",
    "# Get rid of a Matplotlib deprecation warning\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We want to discretize a homogeneous subsurface model with x = 1000 m and z = 500 m and a spatial grid point distance dh = 50 m."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "code_folding": []
   },
   "outputs": [],
   "source": [
    "# Define FD grid parameters\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define no. of grid points in x- and z-direction and grid point coordinates\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nice, next we create and visualize a homogeneous Vp-model discretized on this Cartesian grid."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Hmm, looks a little bit boring, so let's build a more interesting model - e.g. a layer over a half-space:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "... or an even more complex layer interface:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There seem to be something wrong with this model. While there is an undulating layer interface, it does not correctly follow the defined sine-shape. The reason is a spatial undersampling of the model. Let's try a coarser model discretization:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks better, but the interface is still quite blocky. This is a serious disadvantage of  the Cartesian grid and therefore also for the finite-difference approach. We get these \"stair-case\" discretization errors, which become a problem in seismic modelling, because they lead to artificial diffractions. Try to minimize the stair-case artifacts by using an even finer model discretization."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## We learned:\n",
    "\n",
    "* Creating 2D Cartesian grids and P-wave velocity models\n",
    "* Visualization of discrete models\n",
    "* The Cartesian grid has serious limitations to accurately discretize complex layer interfaces which not coincide with the Cartesian coordinate axes"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
